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Probing the lowest energy configuration of a complex system by quantum anneal- 
ing was recently found to be more effective than its classical, thermal counterpart. 
Comparing classical and quantum Monte Carlo annealing protocols on the random 
two-dimensional Ising model we confirm the superiority of quantum annealing relative 
to classical annealing. We also propose a theory of quantum annealing, based on a 
cascade of Landau-Zener tunneling events. For both classical and quantum annealing, 
the residual energy after annealing is inversely proportional to a power of the logarithm 
of the annealing time, but the quantum case has a larger power which makes it faster 

The annealing of disordered and complex systems towards their optimal (or lowest energy) state 
is a central problem in statistical physics, with impact in a large variety of areas. The unknown 
ground state of a system can be approximated by slow-rate cooling of a real or fictitious temperature: 
the slower the cooling, the closer the approximation. Although this kind of standard classical 
annealing (CA) has been extensively investigated over the last two decades and is routinely used 

in a variety of technological applications, such as chip circuitry design, the premium on any alternative 
better optimization algorithms would certainly be enormous. 

Recent results of Brooke et al. on the spin 1/2 disordered Ising ferromagnet LiHoo.44Yo.56F4 

suggested however that a different, quantum annealing (QA) procedure works surprisingly better than 
classical annealing. In QA, temperature is replaced by a quantum mechanical kinetic Hamiltonian 
term - in the specific case a transverse magnetic field T mixing the up and down spin states at each 
site. Initially the quantum perturbation starts out so large in magnitude as to completely disorder 
the system even at zero temperature. When the transverse field is subsequently reduced to zero at 
some slow rate 1/t, the system is "annealed" towards its ground state, much in the same way as 
when its temperature is reduced to zero in CA. The question is which of the two, CA or QA, works 
better, and how and why. Experimental comparison of the properties displayed by the same system 
transported from the same initial state A - a classical high-T state - to the same nominal final state 
B - a low-T glassy state - through two different routes in the [T,r] plane, presents evidence that QA, 
the "quantum route" from A to B, yields with the same "cooling" rate, a state B apparently closer 
to the ground state than CA, the classical one. The data however do not clarify how, and even less 
why, that should be so. 

Theoretical suggestions and exemplifications of QA made by various groups over the past decade 
P~|l0| have stimulated considerable interest in understanding the mechanism of QA better. A theo- 
retical discussion of the relative merits of CA and QA is therefore desirable. For this, it is necessary to 
carry out a direct comparative test on a sufficiently representative benchmark system, such as a spin 
glass, and to lay the bases of a theory of the processes underlying QA. The issues are pressing, both 
because the physical underpinnings of QA call to be explored, and because of the practical potential 
of QA in the fields of optimization in complex systems, should QA turn out to be (as recently shown 
in a protein folding model ||) actually superior. Our work is meant as a step aimed at filling these 
gaps. 

En route, open issues are found even in the context of plain CA, where the very rate of decay 
of the residual energy above the actual ground state energy as a function of the annealing rate 
1/t is controversial. Whereas general theoretical arguments by Huse and Fisher H predict a slow 
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logarithmic convergence e res (r) = Efi na i(r) — E G g ~ log ? (t), with C, < 2, early simulations fLl| |, but 



also more recent studies 12 Jl3], favor a different form, such as power-law, e res (r) ~ r or stretched 



exponential. The question remains whether the discrepancy between simulations and theory is real, 
or only apparent. 

Our work proceeded in three steps. First, we chose a benchmark system, the Ising spin glass, where 
we carried out CA and QA, and compared the results to find that QA is indeed faster. Second, we 
focused on the residual energy in CA, to find deviations from power law decay versus rate 1/t that 
are quite compatible at very slow rates with the Husc-Fisher asymptotically logarithmic decay. Third, 
we built a theory of QA of a spin glass based on the idea of a cascade of level crossings, each with its 
associated Landau-Zener probability to miss the ground state. That theory suggests an asymptotic 
decay of residual energy with QA rate that is again logarithmic as in CA, but governed by somewhat 
different exponents that makes it faster. 

Step 1. Benchmarking quantum versus classical annealing. At the outset, inspired by Brooke 
et al.'s experimental system [|]||, we selected the two-dimensional (2D) random Ising model as an 
appropriate realistic test case. This choice is dictated by the fact the 2D random Ising model is 
technically a polynomial problem Jl4j - where Eqs can be calculated up to sufficiently large lattice 
sizes [jl5|| (thus avoiding an extra fitting parameter in e res (r)) — which is nonetheless of prohibitively 
large complexity for any physical dynamics as a true glass [B 16 1. 



The Edwards- Anderson Hamiltonian of an Ising spin glass in transverse field 

= -£j^f-r$>?, (i) 

(ij) 1 



where nearest-neighbor spins (ij) of a rf-dimensional cubic lattice interact with a random exchange 
coupling Jij, T is the transverse field inducing transitions between the two states, f and |, of each 
spin, and erf, erf are Pauli matrices of the spin 1/2 on site i. The problem is to anneal this system 
as close as possible to its classical, T = 0, ground state. In CA (ljgl, there is no tranverse field and 
no quantum mechanics (r = 0): one starts with a sufficiently high temperature To, which is then 
reduced linearly to zero in a time r. In QA, T is instead fixed to zero or some small value, and one 
starts with a transverse field Tq sufficiently large to throw the system in a "disordered" quantum 
paramagnetic state, decreasing T linearly to zero, again in a time r. Because real time annealing is 
computationally out of the question for the large systems addressed here, we carried out annealing 
as a function, as customary, of the fictitious "time' represented by the number of Monte Carlo steps. 
Our implementation of CA was a standard Metropolis Monte Carlo (MC). That for QA was a path- 
integral Monte Carlo (PIMC) jl7| scheme for a quantum system at a small finite temperature T. The 
2D quantum Ising model is first mapped on a (2+l)D classical model consisting of P copies (Trotter 
replicas) of the original lattice, with a nearest-neighbors uniform ferromagnetic coupling in the third 
(Trotter) direction J x = -(PT/2) logtanh (T/PT), at temperature PT Q. At the beginning of 
the annealing, when T is large, the replicas are only weakly coupled; as T decreases to zero the 
ferromagnetic coupling J 1 - increases, eventually forcing all replicas into the same configuration. At 
the end of either annealing cycle the system, unable to negotiate all barriers in the finite time r, 
remains generally trapped at energy Efi na i — Eqs + £ res, higher than the ground state value Eqs- 
The efficiency of each protocol is measured by the decrease of the average residual energy e res (r) as 
a function of r . 

For a given 2D lattice size L x L, (L up to 80) we took a realization of the random couplings Jy, 
drawn from a flat distribution in the interval (—2,2), and for that we got at the outset the exact 
classical ground state energy Eqs by the Branch and Cut algorithm jljj. Keeping the couplings fixed, 
we then carried out a sufficient number of repeated annealings, (45 for the 80 x 80 lattice), both CA 
and QA. The annealing parameters T (CA) or T (QA) were decreased stepwise from the initial value 
of To — 3 or To = 2.5 down to zero, with a total of r MC steps per spin. In QA we used fixed values 
of PT = 1, 1.5, 2 at several P values, and prepared the initial state (same for all replicas) by classical 
annealing from a temperature of 3.0 down to the corresponding value of PT. In all cases the residual 
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energy e res (r) was calculated by subtracting Eqq from the averaged final annealed energies. 

Fig. 1 shows the residual energy, for both CA and QA for the 80 x 80 lattice, plotted against 
the inverse annealing rate t, actually the actual Monte Carlo computer time. QA appears definitely 
superior to CA, with a lower residual energy for large r. This theoretical finding goes very much in the 
same direction as the experimental evidence of a significantly faster frequency-dependent relaxations 
observed after QA of the disordered magnet Q . The r dependence of our QA data does depend on 
the chosen values of P and T, particularly upon the value of PT, whose optimal value appears to be 
around PT — 1. An increase of P for a fixed value of PT, see inset in Fig. 1, ceases to be effective 
beyond a certain characteristic length (which depends on PT) in the imaginary time direction, . The 
computational cost increases linearly with P, and the choice P = 20 (corresponding to T = 0.05), 
was found to be optimal up to the largest values of r used. Another property (not shown in Fig. 
1) of the CA results is that residual energies obtained for different sizes 32 < L < 80, or even for 
different realizations of the couplings Jy are remarkably size independent and self-averaging, and all 
fall essentially on top of the CA curve in Fig. 1. 

Step 2. Asymptotic behavior of classical annealing. A feature evident in our e res (T) CA data, is 
its gentle but consistent deviation from a pure power law, suggesting serious reconsideration of all 
the earlier power law claims [pH[l2"[ , Since the slope (or apparent power) systematically declines for 
increasing r, it is natural to ask whether it will asymptotically extrapolate to zero in accordance with 
the Huse-Fisher logarithmic law ||. Writing that in the form treJ^ — ^41og(7r) and replacing time 
with number of Monte Carlo steps, we can plot the CA data as in Fig. 2. The extrapolated behavior 
is indeed compatible with a Huse-Fisher straight line. However, as Fig. 2 shows, it proves impossible 
to extract a value for the exponent C, in particular to establish if £ < 2 || is any better, as one could 
have expected. 

Step 3. Landau-Zener theory of quantum annealing. In order to shed some light on the actual 
asymptotic form of residual energy in QA, and eventually rationalize why that might be superior, we 
start off with a cartoon of the instantaneous energy spectrum of (Q) versus V in Fig. 3, suggested 
by small-systems exact diagonalizations. For sufficiently large initial T >> \J%j\ the ground state, 
generally nondegenerate [ fDjfl , must have a finite excitation gap. Imagine following the Schrodinger 
evolution of an initial ground state wavefunction \^>r {t — 0)) while reducing T gradually to zero as 
a function of time [ fL0[ . The instantaneous gap of our disordered magnet will close as T decreases 
through the quantum phase transition at r c p~9] |2l|] . After that, ground state level crossings begin. 
The arrows in the cartoon point to two crossings [really avoided crossings |l8| , the problem possessing 
no symmetry]. Each instantaneous ground state crossing is associated with tunneling of the whole 
system between two valleys - say from a broader but shallower valley to a narrower but deeper one, 
taking place when kinetic energy diminishes - and represents a major crisis in the otherwise quasi- 
adiabatic evolution caused by the time-dependent decrease of T(t). 

For sufficiently slow annealing, each tunneling event can be treated as a Landau-Zener (LZ) problem 
[^2|,^3|, see inset in Fig. 3. The probability P{t) that the system, starting in the lower state \b) at 
high r will continue nonadiabatically onto the higher branch as T is reduced with time is given by 
P(t) = cxp (— t/t c ) where r c , the characteristic tunneling time, is r c = (?iaTo)/(27rA 2 ). Here A is the 
tunneling amplitude between the two states \a) and \b) (whose splitting at crossing is 2| A|), and a is the 
relative slope of the two crossing branches as a function of T |f22| , p3f . One can estimate A ~ e ~ dab ^ ^ , 
where d a b is a suitable distance between states a and b (in the Ising case, the number of spins that are 
flipped in the tunneling process, Nfu p |20j), and £(r) is a typical wavefunction localization length, 
which must vanish as T — > 0, £(r) ~ r^with some exponent </> > 0. The tunneling time becomes 
exponentially large for small T, rr ~ e 2dat, / r *, and an exceedingly small width ~ A of each tunneling 
event justifies treating the multiple crossings as a cascade of independent LZ events. Once the system 
fails, with a probability Pt(t) = e~ r / Tr , to follow the ground state at the LZ crossing occurring at T, 
it will eventually attain an average excitation energy E ex (T) > 0. Letting Z(T)dT be the number of 
LZ crossings which take place between T and T + dT, the average residual energy can be estimated as 
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e res (r) = / dT Z(T)E ex {T) e- T ^ , (2) 
Jo 

where T c marks the first level crossing. The large r behavior of this expression is dominated by the 
r -> behavior of Z(V)E ex (T), and r r . If we assume that, for small T, Z(T)E ex (T) ~ F", and 
Tr ~ e A l v * we get finally a residual energy which vanishes as the inverse power of the rate logarithm 
Crests) ~ log~'* <5A (r), with an exponent Qqa = (1 + ^>)/4>- The exponents ui and cj> are not obvious. 
A semiclassical (WKB) expression for the decay of a wavefunction inside a barrier suggests cf> = 1/2. 
The average excitation energy attained by missing the ground state "track" at P should scale as T 2 
for small T, because all eigenvalues start out as T 2 for T — > 0. The total number of LZ crossings 
occurring from to T should not be larger than the total number of classical states in the energy 
window (Eqs, Eqs + r), which is approximately equal to p(0)T [where p(0) s the density of classical 
states at the ground state energy [jig], so that the density of crossings Z(T — > 0) — > p(0), at most. 
This yields uj = 2 as our most reasonable estimate. 

We conclude that (qa = (1 +u>)/<p can be as large as 6 for a spin glass, and in any case above 
the classical Huse-Fisher bound £ < 2 (3|. Hence, quantum annealing of the Ising spin glass is 
predicted to be again logarithmically accurate, not fundamentally different in that from classical 
annealing. We therefore expect that a quantum computation based on QA will not transform a 
hard nonpolynomial (NP-complete) computational problem into a polynomial one. On the contrary, 
the above reasoning suggests a logarithmically slow annealing to apply also to the present 2D Ising 
case, which is not NP-complete The slowing down effect of the LZ cascade illustrated above 

is particularly severe in problems, like the Ising spin glass we have considered, where the classical 
spectrum has a gapless continuum of excitations above the ground state. Satisfiability problems, for 



which much more encouraging results were recently presented 10 differ from the Ising spin glass in 
that they possess a discrete classical spectrum and a finite excitation gap. We observe that in general 
a gap will cut off the LZ cascade precisely in the dangerous low-F region, and that may eliminate the 
logarithmic slowing down of QA. Nonetheless, even in the gapless case, the advantage of QA over CA 
is far from negligible, due to the generally larger exponent (qa of the logarithm. To get an idea of 
the order of magnitudes involved, consider the relative increase of annealing time (r'/r) needed to 
improve the accuracy of a certain annealing, say with r ~ 10 6 (in appropriate units) by a factor 10. 

In CA (£ = 2), this would require 

( r '/ r ) „ r 10 ' - 1 ~ 10 13 . In QA (f = 6), the same result would 
be accomplished with (r'/r) ~ 10 2 8 , an enormous saving of computer effort. Moreover, the PIMC 
version of QA is easy to implement on a parallel computer, and that provides an extra advantage. 

In summary, our test of QA in the disordered Ising magnet indicates a faster convergence than 
CA, and a time-dependent cascade of Landau-Zener tunneling events across barriers is pinpointed as 
the crucial ingredient of QA. Optimization by QA of a vast variety of problems beyond statistical 
mechanics, of course after a suitable fictitious kinetic energy operator is identified case by case, is an 
open avenue, and stands as a worthy challenge for the future. 
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ful to G. Aeppli, L. Arrachea, J. Berg, C. Micheletti, M. Parrinello, F. Ricci Tersenghi, R. Zecchina, 
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FIG. 1. Comparison of the residual energy per site for an 80 x 80 disordered 2D Ising model after classical 
and quantum annealing. The QA data shown correspond to the optimal value of PT = 1, with T = 0.05 and 
P = 20 Trotter replicas. For fair comparison, the actual inverse annealing rate r used in the QA has been 
rescaled (multiplied by P) so that points at the same r require the same computer time (MCS, Monte Carlo 
steps). The lower residual energy signifies that QA is superior to CA. Inset: r-unrescaled QA data for the 
same system for increasing values of P. Note the satisfactory convergence for P — 20. 
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FIG. 2. The same CA data as in Fig. 1 re-plotted (see text) so as to fall on a straight line if obeying the 
Huse-Fisher logarithmic law. Although the Huse-Fisher form is seen to be asymptotically compatible with the 
data, extraction of a value for the exponent £ is impossible. 
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FIG. 3. Cartoon of the lowest instantaneous eigenvalues of a (finite-size) Ising gl function of the 

transverse field F, or of a generic complex system as a function of its zero-point kinetic energy F. Note the 
two avoided crossing of the ground state, marked by arrows and enlarged in the upper insets. Lower inset: 
Schematic of a Landau-Zener crossing. At each crossing the system will follow adiabatically the ground state 
only if T is reduced sufficiently slowly. The infinite system will exhibit an infinite cascade of crossings as 
T -> 0. 



8 



